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The  physics-based  single  particle  (SP)  model  was  compared  to  the  semi-empirical  equivalent  circuit 
analog  (ECA)  model  to  predict  the  cell  voltage  under  constant  current  charge  and  discharge  for  different 
sets  of  Li-ion  cell  data.  The  parameters  of  the  models  were  estimated  for  each  set  of  data  using  nonlinear 
least  squares  regression.  In  order  to  enhance  the  probability  of  finding  the  global  optima,  a  combination 
of  the  trust  region  method  with  a  genetic  algorithm  was  applied  to  minimize  the  objective  function  (the 
sum  of  squared  residuals).  Several  statistical  quantities  such  as  sum  of  the  squared  errors,  adjusted  R2, 
root  mean  squared  error,  confidence  intervals  of  the  parameters,  and  prediction  bounds  were  included 
to  compare  the  models.  A  significance  test  (t  test)  on  the  parameters  and  the  analysis  of  the  variances  (F 
and  x2  tests)  were  also  performed  to  discriminate  between  the  goodness  of  the  fit  obtained  from  the  two 
models.  The  statistical  results  indicate  that  the  SP  model  superiorly  predicts  all  sets  of  data  compared  to 
the  ECA  model,  while  the  computation  times  of  both  models  are  on  the  same  order  of  magnitude. 
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1.  Introduction 

Today  several  models  for  Li-ion  cells  are  available.  Most  of  these 
models  are  generally  categorized  into  two  main  groups,  physics- 
based  and  empirical.  Physics-based  models  consider  the  physical 
phenomena  which  occur  in  batteries  including  the  material  and 
charge  transfer  processes,  ionic  conduction,  solid  phase  diffusion, 
etc.  Most  of  these  models  are  sophisticated  in  nature  and  time  con¬ 
suming  to  solve  [1-5].  However,  some  efforts  have  been  made  to 
simplify  the  complexity  of  the  physics-based  models  and  reduce 
the  computation  time.  For  example  the  active  material  concentra¬ 
tion  within  the  solid  phase  can  be  approximated  by  the  second 
degree  polynomial  [6-8].  Another  simplification  can  be  obtained 
by  using  the  SP  model  that  represents  each  electrode  as  a  single 
spherical  particle  and  ignores  the  solution  phase  concentration  and 
potential  [9,10].  S.  Santhanagopalan  et  al.  compared  the  two  differ¬ 
ent  physics  based  models;  the  SP  and  the  rigorous  porous  electrode 
models  [11].  The  parameters  of  each  model  were  estimated  and  the 
statistical  F-test  was  performed  to  discriminate  between  the  mod¬ 
els.  A  review  of  the  physics-based  models  comprising  the  rigorous 
pseudo  two  dimensional  (P2D),  the  SP  and  the  porous  electrode 
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models  with  the  polynomial  approximation,  was  presented  by  S. 
Santhanagopalan  et  al.  [12].  The  second  group  of  the  models  for  Li- 
ion  cell  simulation  is  empirical  or  semi-empirical  at  best  [13-15]. 
The  majority  of  these  models  rely  on  the  capacitor/resistor  net¬ 
works  that  are  generally  simple  and  describe  the  physical  processes 
as  being  capacitive  or  resistive  in  nature.  The  original  simplified  cir¬ 
cuit  was  adopted  from  Ong  and  Newman  [16]  by  incorporating  the 
double-layer  capacitance  into  the  model.  M.W.  Verbrugge  and  R.S. 
Conell  implemented  an  equivalent  circuit  model  for  electrochemi¬ 
cal  and  thermal  characterization  of  a  nickel  metal  hydride  traction 
battery  [17].  Model  parameters  were  estimated  by  comparing  the 
average  cell  voltage  and  the  temperature  for  different  charge  and 
discharge  rates.  There  have  been  some  works  that  applied  the  cir¬ 
cuit  analog  approach  for  lead  acid,  nickel  metal  hydride,  Li-ion 
batteries  and  activated-carbon  capacitors  in  adaptive  manner  to 
predict  the  state  of  charge  (SOC),  state  of  health  (SOH)  and  power 
capability  [18-20].  L.  Gao  et  al.  developed  a  model  for  Li-ion  bat¬ 
teries  that  relies  on  an  equivalent  circuit  accounting  for  rate  and 
temperature  dependence  of  the  capacity,  thermal  dependence  of 
the  equilibrium  potential  and  transient  response  [21  ].  Although  use 
of  manufacturers’  data  allows  the  ECA  model  to  have  good  accuracy, 
it  deviates  from  the  experimental  data  at  low  temperature  and  at 
high  discharge  rates.  A  realistic  simulation  of  a  lithium  battery  pack 
was  performed  by  M.  Dubarry  et  al.  [22]  using  an  equivalent  cir¬ 
cuit  approach  for  modeling  of  each  single  cell’s  performance.  They 
showed  that  an  accurate  battery  pack  simulation  can  be  achieved, 
if  the  cell-to-cell  variations  were  taken  into  account.  More  details 
about  the  ECA  model  can  be  found  in  A.J.  Bard  and  L.R.  Faulkner 
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Nomenclature 

an  diagonal  elements  of  the  matrix  A  defined  in  Eq.  (25) 

AR2  adjusted  ft-square 

C  capacitance  [F] 

Cl  confidence  intervals  of  the  parameters 
ce  electrolyte  concentration  [mol  m-3  ] 

ci,max  maximum  solid  phase  concentration  for  each  elec¬ 
trode  (i=p,n)  [molm-3] 
dof  degree  of  freedom 

Ds  i  solid  phase  diffusion  coefficient  of  Li+  or  each  elec¬ 
trode^' =p,n)  [m2  s-1] 

/  the  inverse  of  the  F  cumulative  distribution  function 

F  Faraday’s  constant  [C  mol-1  ] 

H0  null  hypothesis 

Ha  alternative  hypothesis 

lapp  applied  current  [C  s-1  ] 

J  Jacobian  matrix 

Ji  current  density  at  the  surface  of  the  spherical  parti¬ 

cle  for  each  particle  in  electrode  (i=p,n)  [Am-2] 
ki  rate  constant  for  each  electrode  (i=p,n) 

[m2-5  mol-0-5  s-1] 

M  a  large  number 

NCI  normalized  confidence  intervals  of  the  parameters 

rij  number  of  data  points  for  each  dependent  variable 

j 

nm  number  of  points  of  independent  variable 

rip  number  of  parameters 

rip  total  number  of  data  points 
Q.  columbic  capacity  of  the  cell  [Ah] 

PBs  prediction  bounds  for  the  fitted  curves 

R  cell  resistance  [C2] 

R2  ft-squared 

Rct  interfacial  charge  transfer  resistance  [C2] 

Rg  gas  constant  [J  mol-1  K-1  ] 

Ri  radius  particle  for  each  electrode  (i =p,n)  [m] 

r2j-  elements  of  the  correlation  coefficient  matrix 

RMSE  root  mean  squared  error 

s2  approximated  value  for  the  variance 

Sj  electroactive  surface  area  for  each  electrode  ( i = p,n) 

[m2] 

SOQ  state  of  charge  for  each  electrode  (i =p,n) 

SOCceu  cell  state  of  charge 

SOC^u  initial  state  of  charge  for  charge 
SOC^e[[  initial  state  of  charge  for  discharge 
SPB  sum  of  the  absolute  residuals  between  the  predic¬ 
tion  bound  and  the  fitted  model 
SSE  sum  of  squared  errors 

SST  total  sum  of  squares 

t  time  [s] 

T  temperature  [K] 

Uf  open  circuit  potentials  for  each  electrode  (i =p,n)  [V] 
V(bi)  variance  of  parameter  b[ 

Vcen  cell  voltage  [V] 

Vj  electrode  volume  [m3] 

V0  open  circuit  potentials  for  the  cell  [V] 
v  number  of  dependent  variables 

x  predictor  value  defined  in  Eq.  (3 1 ) 

xi,avg  rati°  of  the  solid  average  concentration  to  the  maxi¬ 

mum  solid  concentration  for  each  electrode  (i =p,n) 
xi,surf  ratio  of  the  solid  surface  concentration  to  the  maxi¬ 
mum  solid  concentration  for  each  electrode  (i =p,n) 
rfavgch  ir^tiai  dimensionless  average  concentration  for 
charge  for  each  electrode  (i =p,n) 


v° 

i,avg,  dch 

initial  dimensionless  average  concentration  for  dis¬ 

charge  for  each 

y 

model  prediction 

y 

observed  data 

y 

mean  value  of  the  fitted  responses 

^a,i 

cathodic  transfer  coefficient 

ac,i 

anodic  transfer  coefficient 

Pi 

true  value  of  the  parameter 

m 

overpotentials  for  the  lithium  ion  intercalation  reac¬ 
tion  for  each  electrode  (i =p,n)  [V] 

</>; 

potential  reaction  for  each  electrode  (i=p,n)  [V] 

active  material  volume  fraction  in  the  electrode 

(i=P.n) 

r 

gamma  function 

[23]  and  M.W.  Verbrugge  [24]  books.  None  of  the  above  works 
have  compared  the  physics-based  models  to  the  empirical  or  semi- 
empirical  models  to  investigate  which  models  more  accurately 
predict  experimental  data.  In  this  paper,  a  comparison  of  the  SP 
model  [12]  as  a  physics  based  approach  and  a  simple  ECA  model 
[15]  as  a  semi-empirical  model  was  implemented  by  predicting 
the  charge  and  discharge  voltages  at  different  rates  for  three  sets  of 
data.  The  nonlinear  least  squares  method  was  applied  to  estimate 
the  parameters  of  the  models  by  using  the  Matlab®  Curve  Fitting 
and  Global  Optimization  toolboxes  [23].  Several  statistical  quanti¬ 
ties  and  tests  are  presented  in  order  to  discriminate  between  the 
SP  and  ECA  models. 

2.  Single  particle  model 

A  schematic  of  a  lithium-ion  cell  represented  by  the  SP  model  is 
depicted  in  Fig.  1  with  an  individual  spherical  particle  representing 
each  electrode  as  the  current  being  passed  through  the  electrode  is 
distributed  uniformly  over  all  of  the  particles  in  the  electrode.  The 
main  assumptions  of  the  model  are  as  follows  [12]: 

(1 )  The  concentration  of  the  electrolyte  is  constant  and  uniform  for 
all  time  across  the  cell  sandwich  (cathode,  separator,  anode). 

(2)  The  potential  in  the  solution  phase  is  constant  and  uniform  for 
all  time  across  the  cell. 

(3)  Positive  and  negative  electrode  potentials  depend  on  time  only. 

(4)  There  is  no  side  reaction  in  the  cell. 


Fig.  1.  A  schematic  of  a  lithium-ion  cell  represented  by  the  SP  model  [11]. 
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By  applying  a  volume  average  technique,  Fields  second  law  (par¬ 
tial  differential  equation),  is  simplified  to  predict  the  lithium  ions 
concentration  in  a  spherical  particle.  This  form  contains  an  ordi¬ 
nary  differential  equation  and  algebraic  constraint  equation  (DAE) 
with  the  variables  xi  avg  and  xi<surf  (i=p,n)  for  each  of  the  electrodes 
together  with  their  initial  conditions  as  follows: 


dx; 


i,avg 


dt 


xi,surf  ' 


-3J, 

TRi^i, ms 


-JiRi 


a i,avg  ■ 


5FDs,iCi,  n 


Ic 


Jl  Si 

where  xiAVg  is  the  ratio  of  the  solid  average  concentration  to  the 
maximum  solid  concentration  for  each  electrode  ( ci<max ),  xitSurf  is 
the  ratio  of  the  solid  surface  concentration  to  the  maximum  solid 
concentration  which  is  equal  to  the  state  of  charge  (SOQ)  for  each 
electrode. 

The  Butler-Volmer  kinetic  expression  is  used  to  predict  the  rates 
of  the  lithium  ion  deintercalation  and  intercalation  reactions  for 
each  electrode: 


(1) 


(2) 


(3) 


h 


p  —  Fi(Ci  max  xi, surf  G, max)  (xi,surf^i, max) 


\0.5  0.5 


exp 


VafF 

RgT 


-  exp 


m 


(4) 


The  values  of  the  other  parameters  of  the  SP  model  for  different 
cells  are  given  in  Table  1. 

3.  Equivalent  circuit  analog  model 

Fig.  2  shows  the  schematic  for  an  equivalent  circuit  that  repre¬ 
sents  the  Li-ion  cell  [24]. 

The  single  equation  which  relates  the  cell  current  and  the  volt¬ 
age  takes  the  following  form: 

d(Vceii-Vo)  ,  1 


d  dlapp 

R^r  +  c 


^ —  1  Japp  — 


dt 


+  ^  (Veen -Vo)  (12) 


The  detailed  derivation  is  provided  in  [17].  Since  the  charge  and 
discharge  data  were  obtained  at  constant  current,  the  differential 
equation  for  the  voltage  can  be  written  as: 


d(Vcell-V0)  1  (,  ,  R  \r  1 
dt  -  C  l1  +  i&J  ‘app  -  RaC 


( Vcell  ~  Vo) 


(13) 


where  V0  (open  circuit  potential)  is  a  function  of  the  state  of  charge 
for  the  cell  (SOCceu). 

The  solution  to  the  first  order  differential  equation  (Eq.  (13))  was 
given  in  [17]: 

VCeii  =  V0+IR+Q  exp  +  Japp  Ret  ^1  -  exp  (14) 


The  SOCcen  is  changed  with  time  as  follows: 


SOCcell  —  SOCq  cen  + 


app 


3600 Q 


(15) 


The  overpotentials  for  the  lithium  ion  intercalation  reaction  are: 
m  =  4>i  -  uf  (5) 

The  cell  voltage  is  obtained  by  Eq.  (6): 


VCell  —  0p  -  0n 


(6) 


Since  the  charge  and  discharge  data  were  obtained  at  con¬ 
stant  current,  the  model  does  not  include  the  side  reaction, 
and  the  cathodic  and  anodic  transfer  coefficients  are  the  same 
(aaj  =  aCti  =  o' =  0.5),  the  system  of  differential  algebraic  equations 
(DAEs)  for  the  SP  model  can  be  solved  analytically.  The  analytic 
solution  of  the  differential  equation  for  the  average  concentration 
is  as  follows: 


-3 Ji 


*-i,avg  ■ 


R.Fc, 


-t  +  x" 


(7) 


Surface  concentration  is  obtained: 


xi,surf 


—  Xi 


Jfii 


i,avg 


5  FD<  ,Cj 


(8) 


There  are  two  analytic  solutions  for  the  overpotential  from  the 
Butler-Volmer  equation: 


RgT , 

rii=~E—  In 
Fa 


(Ji  ±  (-4CeF2c2max/<fx2  +  4ceF2C?kfxUurf  +Jf) 


0.5\ 


2 FCq  '^  /q(Cj  ma xX-i,surf  )  (Q, 


;  ci,maxxi,surf) 


,0.5 


(9) 


where  only  the  positive  sign  yields  the  real  value  (physically  mean¬ 
ingful).  Finally,  the  electrode  potentials  can  be  evaluated  by  adding 
the  open  circuit  potentials  to  the  overpotential  obtained  by  Eq.  (9): 

4>i  =  iii  +  u?(xUsurf)  (io) 

The  following  10  parameters  of  the  SP  model  were  selected  to 
fit  the  data: 


where  SOC0)Ce//  is  the  initial  SOC  for  the  cell.  Note  that  the  current 
is  positive  while  charging  the  cell  and  negative  when  the  cell  is 
discharged.  The  open  circuit  potential,  V0,  for  the  cell  is  calculated 
as  the  difference  of  the  open  circuit  potentials  for  the  electrodes: 

Vo  =  Uep(SOCp)-U0n(SOCn)  (16) 

In  order  to  evaluate  the  open  circuit  potentials  for  the  elec¬ 
trodes,  the  state  of  charge  of  the  cell  must  be  correlated  to  the  SOC 
of  each  electrode.  The  values  of  electrode  SOCs  are  known  at  the 
fully  charged  state,  i.e.  SOCceu  =  1.0,  and  the  fully  discharged  state, 
i.e.  SOCcen  =  0-0-  Thus,  the  SOC  for  each  electrode  at  a  particular  cell 
SOC  can  be  determined  by  using  a  linear  interpolation. 

The  following  6  parameters  of  the  ECA  model  were  considered 
to  be  estimated  to  fit  the  data: 

b  =  {C,Q,R,Ra,SOC£cell,SOCt%ll} 

where  SOC^hcell,  SOC^ell  are  the  initial  cell  SOC  for  charge  and  dis¬ 
charge,  respectively. 

4.  Parameter  estimation 

Matlab®  Curve  Fitting  and  Global  Optimization  toolboxes  [25] 
were  used  to  estimate  the  parameters  of  the  SP  and  the  ECA  models 
using  the  nonlinear  least  squares  method.  To  compare  the  models 
statistically,  the  following  are  presented  for  each  model: 

•  Sum  of  squares  due  to  error  (SSE) 

•  R- square  (R2),  adjusted  R- square  (AR2),  root  mean  squared  error 
(RMSE),  confidence  intervals  of  the  Parameters  (C/) 

•  Normalized  confidence  intervals  of  the  parameters  (NCI) 

•  Prediction  bounds  for  the  fitted  curves  (PBs) 

•  Sum  of  the  absolute  residuals  between  the  prediction  bound  and 
the  fitted  model  ( SPB ) 


b  -  {Ds,n,  F)s,p,  Sn,  Sp,  kn,  kp,  xn,avg,ch’ 


x°  x° 

A^n,avg,dch,  p,avg,ch'> 


v° 

p,avg,dch 


} 


(ii) 


Also  the  statistical  t- test  was  performed  to  show  the  significance 
of  the  parameters  of  the  models  and  in  order  to  test  the  adequacy 
of  the  models,  the  F-test  and  the  x2-test  are  presented. 
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Table  1 

SP  model  constant  parameters  in  parameter  estimation  for  different  cells. 


Parameters 

USG 

Quallion 

P2D 

Value 

Unit 

Value 

Unit 

Value 

Unit 

Cp,max 

51,555 

mol  m3 

50,991 

mol  m3 

51,555 

mol  m3 

Cn,max 

30,555 

mol  m3 

30,555 

mol  m3 

30,555 

mol  m3 

Ce 

1000 

mol  m-3 

1000 

mol  m3 

1000 

mol  m3 

RP 

le-5 

m 

5e-6 

m 

8.5e-6 

m 

Rn 

2.5e-5 

m 

1.25e-5 

m 

12.5e-5 

m 

&a,i 

0.5 

- 

0.5 

- 

0.5 

- 

ac,i 

0.5 

- 

0.5 

- 

0.5 

- 

T 

298.15 

K 

298.15 

K 

298.15 

K 

R 

8.3143 

J  mol-1  I<-1 

8.3143 

J  mol-1  K  1 

8.3143 

J  mol-1  K-1 

F 

96,487 

Cmol-1 

96,487 

Cmol-1 

96,487 

Cmol-1 

5.  Nonlinear  least  squares 

Because  the  SP  and  the  ECA  models  are  nonlinear  in  their  param¬ 
eters,  the  nonlinear  least  squares  method  was  applied  to  minimize 
the  SSE  between  the  experiments  and  the  model  given  by  [26]: 


®=EE^i-^)2 

j= 1 


(17) 


where  v  and  rij  are  the  number  of  dependent  variables  and  the  num¬ 
ber  of  data  points  for  each  variable,  respectively;  and  y?.,  are 
the  observed  response  and  the  fitted  response  values,  respectively. 
The  objective  function,  SSE ,  is  minimized  using  the  trust  region 
method  (Curve  Fitting  toolbox)  and  the  genetic  algorithm  (Global 
Optimization  toolbox)  in  the  following  procedure: 

1.  Initial  guess,  lower  and  upper  bounds  for  the  parameters  are 
selected. 


-  Anode  porous  electrode  Separator  Cathode  porous  electrode  + 


zft  'Vf'KS 

>  f‘  J'F  '* 

SCy  ydJr 

J#  Tp  Jr_ 

jmr  *  *  "g/r  9- 

9  .» 

1 

Solid 
«  phase 

i 

i 

i 

i 

i 

i 

electrolyte 

resistance 


electrolyte 

resistance 


Fig.  2.  Equivalent  circuit  for  correlating  the  lithium-ion  cell  behavior  [24]. 
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2.  The  trust  region  method  tries  to  minimize  the  objective  func¬ 
tion  using  the  initial  guess.  Note  that  the  Levenberg-Marquardt 
method  does  not  handle  the  bound  constraints. 

3.  The  resulting  point  from  Step  2  is  considered  as  one  of  the  chil¬ 
dren  of  the  initial  population  in  the  genetic  algorithm  while  the 
other  children  are  generated  randomly  in  the  feasible  domain  for 
the  parameters.  Afterward,  the  genetic  tried  to  find  the  param¬ 
eters  with  less  SSE. 

4.  If  the  objective  function  is  improved  in  Step  3  by  a  defined  tol¬ 
erance  (Simp  =  10-6),  the  optimal  parameters  obtained  in  Step  3 
are  picked  as  the  initial  guess  for  Step  2;  otherwise  terminate 
the  minimization  algorithm  and  print  the  results. 


R2  is  expressed  as  [27]: 


Rz  =  1  - 


SSE 

SST 


where  SST,  the  total  sum  of  squares,  is  defined  as  [27]: 

nT  2 

SST  =  ]T(y*-y) 


1=1 


(18) 


(19) 


where  Cov(bi,  bj ),  V(bt)  are  the  covariance  of  bz  and  bj  and  the  vari¬ 
ance  of  b[  is  obtained  as  follows: 


Cov(bj ,  bj)  =  s2aZj- 
V(bi )  =  s2a;i 


By  substituting  Eq.  (28)  in  Eq.  (27),  the  following  is  obtained: 


(29) 


where  ay  is  the  (ij)  element  of  the  matrix  A  defined  in  Eq.  (25). 
The  higher  the  correlation  between  two  parameters,  the  closer  the 
absolute  value  of  rtj  is  to  1 .0. 

We  defined  the  normalized  parameter  confidence  interval  (NCI) 
in  this  paper  in  order  to  compare  the  confidence  interval  for  differ¬ 
ent  parameters  and  also  for  different  models: 


±ffl-0.05/2)SV/^ 

NCIi  =  — ^ - j2-L -  (30) 

bi 

The  prediction  upper  and  lower  bounds  ( PBUL )  for  the  fitted 
curve  and  for  all  predictor  values  (x)  are  given  by  [27]: 


where  nT  and  y  are  total  number  of  data  points  and  the  mean 
value  of  the  fitted  responses,  respectively.  Note  that  because  all 
the  dependent  variables  are  cell  voltages,  we  define  a  unique  R2  for 
all  the  curves.  Since  the  number  of  parameters  for  the  SP  and  the 
ECA  are  different,  the  AR2  and  the  RMSE  are  also  presented  in  this 
work  as  follows  [27]: 


ar2  =  i  (nr-l)SS£ 

(20) 

(dof -\)SST 

[SSE 

RMSE  =  a  — 

V  dof 

(21) 

dof  =  rip  -np 

(22) 

where  dof  and  rip  are  the  degree  of  freedom  and  the  number  of 
parameters,  respectively. 

The  95%  confidence  interval  for  the  parameter  bz  is  constructed 
as  [26]: 


PBul  =  y  ±fs\JxanXT  (31) 

where/is  the  inverse  of  the  F  cumulative  distribution  function.  The 
PBs  together  with  the  fitted  curves  for  the  SP  and  the  ECA  models 
are  presented.  However,  the  sum  of  the  absolute  residuals  between 
the  lower  bounds  (lower  or  upper  bounds  results  in  no  difference) 
and  the  fitted  model  is  also  calculated  to  compare  the  bounds  of 
the  models  quantitatively: 

V  nj 

SpB  =  J2J2{PBL-yy)  (32) 

j= 1  »=l 

6.  t- test 

A  significant  test  (t- test)  presented  in  [26]  is  to  test  the  null 
hypothesis  that  any  one  of  the  parameters  might  be  equal  to  zero. 
The  null  hypothesis  can  be  stated  as  [26]: 


C/z  -  ±t(i_0.05/2)f 


(23)  H0  :  Pi  =  0 


where  t(i_0.05/2)  *s  a  value  of  Student’s  t-distribution  with  dof 
degrees  of  freedom  and  95%  confidence,  and  s2,  approximated  value 
for  the  variance,  is  given  by  [26]: 


s 


2 


SSE 

dof 


(24) 


azz  in  Eq.  (23)  are  the  diagonal  elements  of  the  following  matrix 
[26]: 


*■(§“) 


(25) 


where  J  is  the  Jacobian  matrix  defined  at  each  point  as  follows: 


and  the  alternative  hypothesis  as 


Ha  :  ^  #  0 


This  is  true  when  the  parameter  appears  in  the  model  as  the 
parameter  vanishes  if  the  value  is  set  to  zero  such  as  polynomial 
coefficients.  Thus,  when  the  test  verifies  the  null  hypothesis,  the 
parameter  can  be  removed  from  the  model.  However,  this  is  not 
the  case  when  the  parameter  is  the  root  of  a  denominator  or  the 
root  of  a  natural  logarithm,  etc.  For  the  former,  /3Z  needs  to  be  a 
large  number  while  for  the  latter  is  one.  For  example,  consider 
the  following  model  with  two  parameters  (bhb2): 

y=  f  +ln(b2) 

b\ 


dSSE 


(26) 


The  elements  of  the  correlation  coefficient  matrix  are  calculated 
by  Eq.  (27)  to  show  how  the  parameters  of  the  models  are  correlated 
[26]: 


Cov{bj,  bj) 


(27) 


The  null  hypothesis  can  be  stated  as 
H0  :  ft  =  M,  =  l 
and  the  alternative  hypothesis  as 
HA  :  M,  /32  /  1 

where  M  is  a  large  number.  As  a  result,  when  the  null  hypothesis  is 
accepted,  the  parameter  can  be  eliminated  in  the  model. 
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Table  2 

Optimization  results  for  fitting  the  USG  and  P2D  data  using  the  SP  and  the  ECA  models. 


8455 


USG  using  the  SP  model 

USG  using  the  ECA  model 

P2D  using  the  SP  model 

P2D  using  the  ECA  model 

Initial  guess/method 

SSE 

Initial  guess/method 

SSE 

Initial  guess/method 

SSE 

Initial  guess/method 

SSE 

Initial  guess 

1060570 

Initial  guess 

295.109 

Initial  guess 

157617 

Initial  guess 

485.126 

Trust  region 

959662 

Trust  region 

2.6536 

Trust  region 

155556 

Trust  region 

0.70651 

Genetic 

0.6043 

Genetic 

2.6536 

Genetic 

0.03437 

Genetic 

0.61433 

Trust  region 

0.4963 

Trust  region 

0.01351 

Trust  region 

0.61368 

Genetic 

0.4963 

Genetic 

0.01351 

Genetic 

0.61368 

The  t  value  is  calculated  for  each  parameter  using  the  following 
equation: 


Fi_Q,(do/i,  do/2)  value  can  be  calculated  by  solving  the  following 
equation: 


bj-Pi 


(33) 


If  this  value  lies  within  the  region  of  acceptance  given  by  the  t 
distribution  at  the  required  confidence  interval  level,  then  the  null 
hypotheses  is  accepted. 


Fi-a(doft,dof. 2) 


/ 


/(t,  do/i ,  dof2)dt  =  1  -  a 


0 


f(t,  dofudof2)m 


F((do/i  +do/2)/2) 
r(do/i/2)r(do/2/2) 


dofi 

2 


(36) 

t(doh-2)l2 

[\+(d0hld0f2)tf0J'+d0J2>/Z 


where  r  is  the  gamma  function. 


7.  F-test 


8.  Chi-squared  (x2)  test 


In  order  to  test  the  adequacy  of  the  fit  of  the  SP  and  the  ECA 
models,  the  F-test  was  performed  [26].  The  total  sum  of  residuals 
is  the  sum  of  the  error  due  to  the  lack  of  fit  and  the  experimental 
error.  The  variances  of  these  errors  are  obtained  using  the  following 
equations: 


S 


2  _ 
1  ~~ 


doh 

nm 

X>*-y*)2 


s2  =  i=! _ 

2  dof2 

do/!  =  nnr  -  np 
dof2  =  nT-  n„r 


(34) 


where  y*  is  the  mean  value  of  each  group  of  repeated  experiments, 
nm  is  the  number  of  points  of  independent  variable  (the  total  non 
repeated  points),  and  nT  is  the  total  number  of  points.  The  ratio  of 
the  variances,  s\\s\,  has  an  F distribution  with  do/i  and  dof2  degrees 
of  freedom.  For  a  good  fit,  the  following  inequality  must  be  true: 

s2 

-1  <F1_0(do/,,do/2)  (35) 

S2 


In  order  to  show  the  goodness  of  the  fit,  the  chi-squared  (x2) 
statistical  test  is  also  accomplished  [28].  The  quantity  known  as 
chi-squared  (x2)  is  calculated  as  follows: 


(y»  —  y*) 


va * 


1=1 

n* 

va *  =  ^(y*  -  y* 

j= 1 


(37) 


This  quantity  has  a  x2  distribution  with  dof\  degrees  of  freedom 
and  to  have  a  good  fit,  the  following  muse  be  satisfied: 

X2<Xdo/l,1_„  (38) 

where  xdo^  i_a  *s  obtained  by  solving  the  following  equation: 


00 


/ 


/(/W  =  l-« 


^ dofi  >  1  ~a 

Kx2)  = 


e-/2ll(yl  j(cf0/l/2)-  l 

2d°fi  /2  r  ^ j 


(39) 


Table  3 

Statistical  results  for  the  USG  data  fitting  using  the  SP  model. 


Parameter 

10  Parameter  model 

8  Parameter  model 

Value 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

Ds  n  (m2  s-1) 

1.3863e-14 

-27691 

Yes 

1.3862e-14 

±4.0258e-15 

±0.29042 

-48678 

Yes 

Ds,p  (m2  s-1) 

9.9984e-ll 

-1.9e-6 

No 

- 

- 

- 

- 

- 

S„(m2) 

0.14477 

- 

- 

0.14489 

±0.021376 

±0.14753 

- 

- 

Sp(m2) 

0.18868 

- 

- 

0.18867 

±0.0018607 

±0.0098621 

- 

- 

kn  (m2-5  mol--5  s-1) 

9.9991e-6 

-2.3e-9 

No 

- 

- 

- 

- 

- 

kp  (m2-5  mol--5  s-1) 

4.7568e-12 

-2.6e7 

Yes 

4.7585e-12 

±1.9004e-13 

±0.039937 

-1.03e8 

Yes 

0.11119 

7.1139 

Yes 

0.11136 

±0.015185 

±0.13636 

14.244 

Yes 

n,avg,cn 

0.83869 

29.497 

Yes 

0.83814 

±0.053812 

±0.064204 

30.491 

Yes 

n,avg,acn 

0.99534 

6794.1 

Yes 

0.99534 

±0.00014166 

±0.00014232 

13757 

Yes 

p,avg,cn 

0.5313 

171.59 

Yes 

0.53128 

±0.004497 

±0.0084645 

231.11 

Yes 

p,avg,acn 
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9.  Results  and  discussion 

9.1.  USG  Gen  1.0  experimental  data 

The  experimental  data  consists  of  the  cell  voltage,  current  and 
the  temperature  for  several  charge/discharge  cycles  of  the  lithium- 
ion  pouch  cell.  The  rated  capacity  of  the  cell  is  341  mAh.  The  active 
materials  of  the  positive  and  negative  electrodes  are  lithium  cobalt 
oxide  (LiCo02)  and  mesocarbon  microbead  (MCMB)  2528,  respec¬ 
tively.  The  three  first  charge  curves  at  0.2  C-rate  and  three  discharge 
curves  at  different  rates,  0.5, 0.686,  and  1 .029  C-rates  were  selected 
to  estimate  the  parameters  of  the  SP  and  the  ECA  models.  The  anode 
and  cathode  open-circuit  potentials  as  functions  of  the  SOCs  are 
expressed  as  [5]: 

U„(SOC„ )  =  .7222  +  .1387 SOC„  +  .029 SOC’/2  -  2171  + 

SULn  SOCn' 

+  .2808e(°-9_15SOCn)  -  .7984e(0-4465SOCn_0-4108) 


Up(SOCp)  = 


-4.656  +  88.669SOC2  -  401.119 SOC4  +  342.909 SOCj] 
-1  +  18.93350Cp  -  79.532SOC4  +  37.31  lSOC^  - 


It  is  assumed  that  the  anode  and  cathode  SOCs  vary  linearly  with 
the  cell  SOC  as  already  mentioned.  Fig.  1  in  [29]  was  used  to  obtain 
the  following  linear  functions: 


f  SOCa  =  0.79SOCcen  +  0.01 
\  SOCc  =  0.97  -  0. 51  SO Cce// 


9.1.1.  SP  model 

The  SP  model  includes  the  following  10  parameters: 


b  -  {Ds,n,  Ds,p,  Sn,  Sp,  kn,  kp ,  X°  avg  ch,  ^  avg  dch,  avg)dl,  xp,avg,dch} 


The  diffusion  term  in  the  model  can  be  removed  when  the  diffu¬ 
sion  coefficient  becomes  large.  As  a  result,  fa  in  the  t  value  equation 
is  set  to  the  upper  bound  (10-5  m2  s-1)  used  in  the  optimization. 
The  active  surface  area  (S*)  is  also  required  to  be  a  large  number 
for  elimination  in  the  model.  However,  the  upper  bound  for  this 
parameter  is  limited  to  the  theoretical  value  that  cannot  be  consid¬ 
ered  as  the  value  for  fa.  Thus,  the  t  value  is  not  defined  for  the  Sp 
In  order  to  eliminate  the  kinetic  parameters  (/q)  in  the  model,  the 
expression  in  the  natural  logarithm  of  the  overpotential  equation 
is  needed  to  be  one.  This  occurred  when  the  /q  limit  to  the  infin¬ 
ity,  which  means  there  is  no  reaction  resistance  for  intercalation 
or  deintercalation  of  Li+  ions  in  the  electrodes.  Thus,  fa  needs  to  be 
a  large  number.  The  upper  bound  (10-1°  m2-5  mol-0,5  s-1 )  used  in 
the  optimization  for  the  /q  is  considered  as  fa.  The  fa  values  for  the 
initial  SOCs  are  zero  as  these  parameters  are  linear  in  the  model. 
Table  2  presents  the  objective  function  improvement  by  using  the 
optimization  algorithm.  The  genetic  algorithm  greatly  improved 
the  objective  function  after  the  trust  region  method  as  shown  in 
Table  2.  The  parameters’  values  and  the  corresponding  t  values 
are  given  in  Table  3.  The  table  indicates  that  the  cathode  diffu¬ 
sion  coefficient  and  the  anode  rate  constant  do  not  affect  the  fitting 
significantly.  Thus,  the  regression  was  repeated  for  the  model  with 
8  remaining  parameters.  In  the  new  model,  since  there  is  no  diffu¬ 
sion  resistance  for  the  positive  electrode,  the  surface  concentration 
is  the  same  as  the  average  concentration  and  the  negative  electrode 
potential  is  equal  to  the  open  circuit  potential  because  there  is  no 
reaction  resistance: 


xp,surf  ~  xp,avg 
0n  =  Un(xn,surf ) 


The  initial  state  of  charge  of  the  positive  electrode  returned 
from  the  parameter  estimation  routine  is  larger  than  expected 
(0.99  <  x9  avg  ch  <  1 .0).  The  large  value  is  likely  due  to  the  fact  that 
the  SP  model  does  not  take  into  account  the  potential  drop  in  the 
electrolyte. 

Table  3  consists  of  the  parameters’  values,  confidence  intervals, 
normalized  confidence  intervals,  and  t  values.  The  same  SSE  for 
the  10  parameter  SP  model  was  obtained.  The  fitted  SP  model  and 
the  prediction  bounds  for  different  charge  and  discharge  rates  are 
shown  in  Fig.  3a  and  b,  respectively.  In  order  to  check  the  validity 
of  the  values  of  the  active  surface  areas  obtained  by  the  parameter 
estimation,  the  active  material  volume  fractions  in  the  electrodes 
are  calculated  as  follows: 

p 

*  “  31/j 

The  projected  areas  for  the  positive  and  negative  electrodes  are 
0.0213  and  0.0226  m2,  respectively.  Since  the  thicknesses  are  53 
and  72  [xm  for  the  positive  and  negative  electrodes  respectively, 


-  462. 471  SOC8  +  433.434SOCp° 
73.083 SOC8  +  95.96 SOC40 


the  active  material  volume  fractions  were  calculated  as  0.558  for 
the  positive  and  0.743  for  the  negative. 

9.1.2.  ECA  model 

The  following  6  parameters  were  used  in  the  ECA  model  to  fit 
the  data: 

{C,<lR,Ra,SOC*aa,SOC™ett} 

The  value  of  fa  for  capacitance  C  is  a  large  number  (e.g.  upper 
bound  of  the  optimization  (1 0,000  F))  because  as  capacitance  C 
limits  to  the  infinity,  the  corresponding  term  in  the  model  can 
be  eliminated.  For  other  parameters,  fa  is  zero.  The  optimization 
results  are  shown  in  Table  2.  In  this  case,  the  genetic  did  not  improve 
the  SSE  after  the  trust  region  method.  Table  4  presents  the  least 
squares  results  including  the  parameters’  values  and  the  t  values. 
The  t  value  for  resistance  R  lies  within  the  region  of  acceptance 
of  the  null  hypothesis  given  by  the  t  distribution  at  the  95%  con¬ 
fidence  interval.  This  means  that  resistance  R  can  be  removed  in 
the  model  and  as  a  result  we  need  to  repeat  the  least  squares  algo¬ 
rithm  for  the  model  with  remaining  5  parameters.  The  results  are 
given  in  Table  4.  It  can  be  seen  that  the  t  value  of  the  interfacial 
charge  transfer  resistance,  Rct,  is  enhanced  more  with  respect  to 
the  other  parameters  when  the  resistance  R  is  eliminated.  Fig.  4a 
and  b  illustrates  the  fitted  ECA  model  and  the  prediction  bounds 
for  the  charge  and  discharge  rates,  respectively. 

9.1.3.  Comparison 

Table  5  presents  the  statistical  quantities  resulting  obtained  by 
fitting  the  USG  data  using  the  SP  and  the  ECA  models.  Based  on  the 
results  in  the  table,  it  is  obvious  that  the  SP  model  predicts  the  USG 
data  much  superior  than  the  ECA  model. 

9.2.  Quallion  Gen  2.0  experimental  data 

One  charge  curve  and  one  discharge  curve  both  at  0.393  C-rate 
(0.9825  A)  were  selected  among  several  charge/discharge  cycles’ 
experimental  data  obtained  from  the  natural  graphite  (Mag-1 0)/Li 
half-cell,  LiNio.8Coo.15Alo.05O2/Li  half-cell.  The  rated  capacity  of  the 
cell  is  2.5  Ah.  The  open  circuit  potentials  versus  the  SOC  of  the  active 
material  in  the  positive  and  negative  electrodes  have  been  reported 
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Fig.  3.  The  fitted  SP  model  with  8  parameters  and  the  prediction  bounds  for  USG  data  (a)  charge  and  (b)  discharge. 


Table  4 

Statistical  results  for  the  USG  data  fitting  using  the  ECA  model. 


Parameter 

6  Parameter  model 

5  Parameter  model 

Value  t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

C  (F) 

1795.6  -16.57 

Yes 

1793.8 

±514.3 

±0.2869 

'  -31.278 

Yes 

Q(Ah) 

0.49344  364.38 

Yes 

0.49344 

±0.0026525 

±0.0053757  364.6 

Yes 

R(Q) 

1.170e-8  3.727e-7 

No 

- 

- 

- 

- 

- 

RctW 

0.27783  8.5218 

Yes 

0.2778 

±0.012425 

±0.044728  43.819 

Yes 

SOCceiicn 

0.037499  43.383 

Yes 

0.037503 

±0.0016396 

±0.043728  44.821 

Yes 

SOCCeiidch 

0.82258  370.92 

Yes 

0.82258 

±0.0043451 

±0.0052822  371.05 

Yes 

in  [30].  The  following  linear  functions  for  the  anode  and  the  cathode 

both  electrodes  and  all  initial  SOCs  do  not  influence  the  fitting 

SOCs  variations  with  the  cell  SOC  obtained  based  on  Fig.  3  in  [30]: 

significantly.  However,  the  correlation  coefficient  matrix  for  the 

parameters  states  that  the  SOCs  for  each  electrode  are  correlated 

/  SOCa  = 

0.7388SOCce/j  +  0.0252 

with  the  electrode  diffusion  coefficient  as  shown  in  Table  7.  As  a 

\  SOCc  = 

— 0.62  +  0.96 SOCceu 

result,  we  repeated  the  regression  with  8  parameter  model  that 

excludes  the  diffusion  coefficients  because  their  t  values  are  less 

than  the  others.  Thus,  the  Li  ion  surface  concentrations  are  assumed 

9.2.1.  SP  model 

to  be  equal  to  the  average  concentration  for  both  electrodes.  The 

Table  6  lists  the  values  of  the  parameters  and  the  t  values 

statistical  results  are  presented  in  Table  6.  Fig.  5a  and  b  shows  the 

obtained  by  fitting  the  Quallion  data  using  the  SP  model  with  10 

fitted  SP  model  and  the  prediction  bounds  for  charge  and  discharge 

parameters.  The  table  indicates  that  the  diffusion  coefficients  for 

rates,  respectively. 

Table  5 

Comparison  of  the  SP  model  with  the  ECA  model  for  fitting  different  data. 

Data 

Model 

#  of  the  parameters 

dof 

SSE 

R2 

AR2 

RMSE 

SPB 

s2 

4  X2  value 

S2 

USG 

SP 

8 

2704 

0.49625 

0.99214 

0.99212 

0.013547 

72.147 

- 

ECA 

5 

non 

2.6536 

0.95798 

0.95792 

0.031309 

166.65 

- 

Quallion 

SP 

8 

100 

0.01414 

0.99911 

0.99905 

0.011891 

2.6431 

- 

ECA 

5 

103 

0.24647 

0.98456 

0.98396 

0.048917 

10.713 

- 

ECAa 

8 

100 

0.039934  0.9975 

0.99732 

0.019983 

4.4325 

P2D 

SP 

9 

231 

0.01351 

0.99864 

0.99859 

0.007647 

3.6814 

1.6582  344.05 

ECA 

5 

235 

0.61368 

0.93808 

0.93702 

0.051102 

24.404 

74.041  14210 

ECA  with  quadratic  Rct. 
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Fig.  4.  The  fitted  ECA  model  with  5  parameters  and  the  prediction  bounds  for  the  USG  data  (a)  charge  and  (b)  discharge. 


Table  6 

Statistical  results  for  the  Quallion  data  fitting  using  the  SP  model. 


Parameter 

10  Parameter  model 

8  Parameter  model 

Value 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

DSji  (m2  s_1 ) 

3.7077e-ll 

-8.1e-9 

No 

- 

- 

- 

- 

- 

Ds,p  (m2  s-1) 

9.999e-ll 

-2e-15 

No 

- 

- 

- 

- 

- 

Sn  (m2) 

0.77424 

- 

- 

0.77423 

±0.0047632 

±0.0061521 

- 

- 

Sp  (m2) 

2.2962 

- 

- 

2.2962 

±0.034462 

±0.015008 

- 

- 

kn  (m2-5  mol--5  s_1) 

1.3572e-10 

-1.8e5 

Yes 

1.3567e-10 

±1.3742e-10 

±1.0129 

-1.4e5 

Yes 

kp  (m2-5  mol--5  s_1) 

1.0191e-12 

-7.2e7 

Yes 

1.0192e-12 

±2.7152e-13 

±0.2664 

-7.3e7 

Yes 

x°  . 

0.045253 

7.2e-6 

No 

0.045282 

±0.0041508 

±0.091666 

21.165 

Yes 

n,avg,cn 

X°  .. 

0.92108 

1.5e-4 

No 

0.92105 

±0.0056974 

±0.0061857 

320.39 

Yes 

n,avg,acn 

X°  . 

0.90801 

1.7e-5 

No 

0.908 

±0.010002 

±0.011015 

179.91 

Yes 

p,avg,cn 

x° 

v.avs.dch 

0.42561 

8e-6 

No 

0.42561 

±0.0090567 

±0.021279 

93.015 

Yes 

Table  7 

Correlation  coefficient  matrix  for  the  Quallion  data  fitting  using  the  SP  model  with  10  parameters. 

Parameters 

Ds,n 

Ds,P 

Sn 

Sp 

kn 

kp 

*i  a 

X2b 

x3c 

*4d 

Ds,n 

1 

0.32 

-5.6e-7 

2.4e-7 

-1.6e-7 

1.42e-7 

1 

-1 

-0.32 

0.32 

Ds,p 

0.32 

1 

-1.3e-7 

3.8e-8 

-1.9e-8 

4.2e-8 

0.32 

-0.32 

-1 

1 

Sn 

-5.6e-7 

-1.3e-7 

1 

-0.4 

0.22 

-0.26 

-4.4e-7 

8.7e-8 

1.6e-7 

-1.6e-7 

sp 

2.4e-7 

3.3e-8 

-0.4 

1 

-0.61 

0.06 

3.9e-8 

-4.5e-8 

-6.4e-8 

5.2e-8 

kn 

-1.6e-7 

-1.9e-8 

0.22 

-0.61 

1 

-0.29 

4.1e-8 

3.1e-8 

3.1e-8 

-1.5e-8 

K 

1.4e-7 

4.2e-8 

-0.26 

0.06 

-0.29 

1 

1.3e-7 

-1.2e-8 

-1.3e-7 

1.2e-7 

*1  3 

1 

0.32 

-4.4e-7 

3.9e-8 

4.1e-8 

1.3e-7 

1 

-1 

-0.32 

0.32 

*2b 

-1 

-0.32 

8.7e-8 

-4.5e-8 

3.1e-8 

-1.2e-8 

-1 

1 

0.32 

-0.32 

*3C 

-0.32 

-1 

1.6e-7 

-6.4e-8 

3.1e-8 

-1.3e-7 

-0.32 

0.32 

1 

-1 

X4d 

0.32 

1 

-1.6e-7 

5.2e-8 

-1.5e-8 

1.2e-7 

0.32 

-0.32 

-1 

1 

b 


d 


X\  =  X°  .  . 

1  n,avg,ch 

*2  =  Xn,ovg,<idi- 

*3 =  Xp,mg,ch ' 

X*=Xp,mg,dC- 
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Fig.  5.  The  fitted  SP  model  with  8  parameters  and  the  prediction  bounds  for  the  Quallion  data  (a)  0.393  C-rate  charge  and  (b)  0.393  C-rate  discharge. 


Table  8 

Statistical  results  for  the  Quallion  data  fitting  using  the  ECA  model. 


Parameter 

6  Parameter  model 

5  Parameter  model 

Value 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly  affecting 
the  fitting? 

C(F) 

128.43 

-53.725 

Yes 

87.275 

±155.58 

±1.7826 

-126.36 

Yes 

Q(Ah) 

2.833 

65.794 

Yes 

2.8369 

±0.078185 

±0.02756 

71.961 

Yes 

R(Q) 

0.031992 

0.92098 

No 

- 

- 

- 

- 

- 

Rct(V) 

0.030438 

0.92165 

No 

0.061486 

±0.014926 

±0.24276 

8.1696 

Yes 

soccellch 

0.081913 

5.5514 

Yes 

0.084271 

±0.026832 

±0.3184 

6.2288 

Yes 

SOECelldch 

0.84882 

65.796 

Yes 

0.84765 

±0.023361 

±0.027559 

71.963 

Yes 

9.2.2.  ECA  model 

The  parameters’  values  with  corresponding  t  values  obtained  by 
fitting  the  Quallion  data  using  the  ECA  model  with  6  parameters  are 
listed  in  Table  8.  Although  the  t  values  for  both  resistances  R  and  Rct 
verify  the  null  hypothesis,  the  correlation  coefficient  matrix  shows 
a  correlation  between  them.  Thus,  the  model  that  excludes  only 
resistance  R  is  used  to  fit  the  data.  The  parameters  values  together 
with  confidence  interval  and  t  values  are  given  in  Table  8.  The  fitted 
ECA  model  and  the  prediction  bounds  for  the  charge  and  discharge 
rates  are  depicted  in  Fig.  6a  and  b,  respectively. 


9.2.3.  Comparison 

The  comparison  between  the  SP  and  the  ECA  models  are  given  in 
Table  5  that  indicates  that  the  SP  model  goodness  of  the  fit  surpasses 
the  fit  for  the  ECA  model. 


9.3.  Pseudo  two  dimensional  (P2D)  model  data 

In  order  to  accomplish  the  F  and  x2  statistical  tests  for  the  SP 
and  the  ECA  models,  a  zero  mean  Gaussian  noise  with  the  speci¬ 
fied  standard  deviation  (e.g.  0.01 )  is  added  to  the  data  obtained  by 
running  the  rigorous  P2D  model  in  COMSOL  3.5a  environment  [31  ] 
under  the  four  different  charge/discharge  rates,  0.1 , 0.25, 0.5  and  1 .0 
C-rate.  Three  replicates  were  generated  at  each  independent  vari¬ 
able  point.  The  same  chemistry  as  the  USG  cell  but  with  different 
dimensions  was  used  and  the  cell  rated  capacity  was  1.656  Ah. 

9.3.1.  SP  model 

Table  9  lists  the  values  of  the  parameters  and  the  t  values 
obtained  by  fitting  the  P2D  data  using  the  SP  model  with  1 0  param¬ 
eters.  The  table  indicates  that  the  diffusion  coefficients  for  the 
positive  electrodes  do  not  affect  the  fitting  significantly.  Thus,  the 


Fig.  6.  The  fitted  ECA  model  with  5  parameters  and  the  prediction  bounds  for  the  Quallion  data  (a)  0.393  C-rate  charge  and  (b)  0.393  C-rate  discharge. 
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Table  9 

Statistical  results  for  the  P2D  data  fitting  using  the  SP  model. 


Parameter  10  Parameter  model 


9  Parameter  model 


Value 

t  value 

Is  parameter 

Value 

Cl 

NCI 

t  value 

Is  parameter 

significantly 

significantly 

affecting  the 

affecting  the 

fitting? 

fitting? 

DSin  (m2  s-1) 

5.5529e-14 

-2.9e5 

Yes 

5.5529e-14(5.5e-14)a 

±6.767e-16 

±0.012186 

-2.9e5 

Yes 

Ds,p  (m2  s-1) 

le-10 

-2e-16 

No 

— (le— 1 1  )a 

- 

- 

- 

- 

Sn  (m2) 

0.68396 

- 

- 

0.68395  (0.78794)a 

±0.010436 

±0.015258 

- 

- 

Sp  (m2) 

1.131 

- 

- 

1.131  (1.1442)a 

±0.010912 

±0.0096481 

- 

- 

kn  (m2-5  mol--5  s-1) 

1.1067e-ll 

-1.4e7 

Yes 

1.1068e-ll  (1.764e-ll)a 

±1.3737e-12 

±0.12412 

-1.4e7 

Yes 

kp  (m2-5  mol--5  s-1) 

4.2239e-ll 

-9.24e5 

Yes 

4.219e-ll  (6.667e-ll)a 

±1.7664e-ll 

±0.41867 

-l.le6 

Yes 

v° 

n,avg,ch 

0.070546 

11.614 

Yes 

0.070486  (0.05)a 

±0.011504 

±0.1632 

11.901 

Yes 

X° 

n,avg,dch 

0.82842 

124.79 

Yes 

0.82845  (0.756)a 

±0.012442 

±0.015018 

131.04 

Yes 

X° 

p,avg,ch 

0.97272 

528.63 

Yes 

0.9727  (0.95)a 

±0.0033922 

±0.0034874 

564.39 

Yes 

X° 

p,avg,dch 

0.46384 

536.91 

Yes 

0.46384  (0.465)a 

±0.0011285 

±0.0024329 

808.09 

Yes 

a  P2D  model  parameters. 


regression  was  repeated  with  9  parameter  model  excluding  the 
positive  electrode’s  diffusion  coefficient.  Table  9  and  Fig.  7a  and  b 
present  the  results  (only  the  data  for  1 .0  C-rate  ( 1 .656  A)  charge  and 
discharge)  obtained  by  fitting  the  P2D  data  using  the  SP  model  with 
9  parameters.  The  genetic  algorithm,  also,  in  this  case  improved 
the  objective  function  significantly  after  the  trust  region  method 
as  shown  in  Table  2.  The  input  parameters  of  the  P2D  model  are 
also  given  in  Table  9  that  indicates  the  parameters  of  the  P2D  and 
the  SP  models  are  in  the  same  order  of  magnitude. 


9.3.2.  ECA  model 

The  parameter  estimation  results  obtained  by  fitting  the  P2D 
data  using  the  ECA  model  with  6  parameters  are  listed  in  Table  10. 
Although  the  t  values  for  both  resistances  R  and  Rct  verify  the 
null  hypothesis,  the  correlation  coefficient  matrix  shows  a  correla¬ 
tion  between  them.  Based  on  the  model  equations,  we  decided  to 
exclude  resistance  R  since  the  resistance  Rct  is  in  the  two  terms  of 
the  model  while  resistance  R  is  only  in  one  term  of  the  equation. 
The  remaining  5  parameters’  values  together  with  their  confidence 
intervals  and  t  values  are  given  in  Table  10.  The  improvement  of 
the  objective  function  is  presented  in  Table  2.  The  fitted  ECA  model 
and  the  prediction  bounds  for  only  1.0  C-rate  charge  and  discharge 
are  depicted  in  Fig.  8a  and  b. 


9.4.  Comparison 

Table  5  compared  the  SP  and  the  ECA  models  for  fitting  the 
P2D  data.  At  the  95%  confidence  interval,  the  values  for  F95%(do/i, 
dof2)  and  x2do^  95%  f°r  the  SP  model  are  calculated  as  1.2009  and 


267.46  respectively  using  the  Matlab®  Statistics  Toolbox  and  for 
the  ECA  model,  the  F95%(do/i,  do/2)  and  xdofA  95%  values  are  271.76 
and  1.1997  respectively.  The  last  two  columns  of  the  table  are  the 
ratio  of  the  model  and  the  experimental  variances  and  x2  values 
for  each  model.  The  data  contradict  the  goodness  of  both  models, 
since  the  ratios  of  the  variances  and  the  x2  values  are  more  than 
the  corresponding  distribution  values.  However,  the  values  for  the 
SP  model  are  much  less  than  the  values  for  the  ECA  model.  The 
other  data  in  the  table  also  verify  the  ascendancy  of  the  goodness 
of  the  fit  for  the  SP  model  with  respect  to  the  ECA  model.  The  com¬ 
putation  times  to  evaluate  the  voltages  at  different  rates  for  both 
models  are  about  0.001  sec  on  a  Dell  Precision  T7500,  with  2  Quad 
Core  2.53  GHz  Zenon  Processors  CPUs  and  12.285  GB  of  RAM. 

10.  ECA  model  with  quadratic  Rct 

In  order  to  compare  the  ECA  model  with  the  same  number  of 
parameters  as  the  SP  model,  it  was  assumed  that  the  interfacial 
charge  transfer  resistance,  Rct ,  is  a  linear  function  of  the  cell  SOC  as 
follows  [17]: 

Rch  =  rch  +  rchS0Ccell  +  rfSOC2cell 
Rdch  =  rdch  +  rfhSOCcell  +  4chSOC2ce[l 

Thus,  the  following  10  parameters  are  required  to  be  estimated: 


{c,  Q,  rf,  rf,  rf 


l  ru 
’  r1 


\  SOC“ 


cell 


SOCdch  \ 
,cellf 


Note  that  the  resistance  R  was  excluded  based  on  the  previous 
results.  The  parameters,  rvalues  and  confidence  intervals  obtained 


Fig.  7.  The  fitted  SP  model  with  9  parameters  and  the  prediction  bounds  for  the  P2D  data  (a)  1.0  C-rate  charge  and  (b)  1.0  C-rate  discharge. 
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Table  10 

Statistical  results  for  the  P2D  data  fitting  using  the  ECA  model. 


Parameter  6  Parameter  model 


5  Parameter  model 


Value 

t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

C  (F) 

21.045 

-2201.4 

Yes 

22.43 

±9.6562 

±0.43049 

-2035.7 

Yes 

Q.(Ah) 

2.0764 

122.54 

Yes 

2.086 

±0.030784 

±0.014757 

133.5 

Yes 

R(Q) 

0.029635 

1.2631 

No 

- 

- 

- 

- 

- 

RctW 

0.010328 

0.45509 

No 

0.038705 

±0.0077621 

±0.20054 

9.8237 

Yes 

SOCceach 

0.91685 

123.19 

Yes 

0.91262 

±0.013383 

±0.014664 

134.35 

Yes 

SOECelldch 

0.037124 

6.1923 

Yes 

0.041068 

±0.010204 

±0.24846 

7.9294 

Yes 

Fig.  8.  The  fitted  ECA  model  with  5  parameters  and  the  prediction  bounds  for  the  P2D  data  (a)  1.0  C-rate  charge  and  (b)  1.0  C-rate  discharge. 


Table  11 

Statistical  results  for  the  Quallion  data  fitting  using  the  ECA  model  with  quadratic  Rct . 


Parameter  10  Parameter  model 


8  Parameter  model 


Value 

t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

Value 

Cl 

NCI 

t  value 

Is  parameter 
significantly 
affecting  the 
fitting? 

C  (F) 

75.467 

-324.51 

Yes 

48.68 

±38.864 

±0.79835 

-508 

Yes 

Q.(Ah) 

3.0363 

191.01 

Yes 

3.0462 

±0.031447 

±0.010324 

192.18 

Yes 

r?  («) 

0.10065 

3.2704 

Yes 

0.1078 

±0.035978 

±0.33376 

5.9443 

Yes 

rch 

' 1 

0.0010993 

0.012323 

No 

- 

- 

- 

- 

- 

rch 

*2 

-0.07627 

-0.77127 

No 

- 

- 

- 

- 

- 

rQCh  W 

0.21198 

31.033 

Yes 

0.20728 

±0.01439 

±0.069423 

28.578 

Yes 

rdch 

1 1 

-8.6635 

-13.777 

Yes 

-7.9307 

±1.2142 

±0.1531 

-12.959 

Yes 

*2 

10.594 

13.18 

Yes 

9.7004 

±1.678 

±0.17298 

11.469 

Yes 

SOCcellch 

0.095094 

5.0206 

Yes 

0.069774 

±0.035081 

±0.33376 

3.946 

Yes 

SOECelldch 

0.7955 

191.01 

Yes 

0.79291 

±0.0081856 

±0.069423 

192.18 

Yes 

for  Quallion  data  fitting  are  presented  in  Table  1 1  that  indicates 
that  the  linear  and  quadratic  coefficients  for  the  charge  Rct  does 
not  affect  the  fitting  significantly.  Thus,  the  charge  Rct  was  kept 
constant  for  the  subsequent  estimation.  The  Rct  linear  coefficient  at 
discharge  (r^ch)  makes  the  resistance  take  negative  values  at  some 
cell  SOC.  In  this  condition  a  small  positive  value  (e.g.  le-6)  is  used 
for  the  resistance.  Other  statistical  quantities  are  given  in  Table  5 
that  shows  the  improvement  of  the  ECA  model.  However,  the  SP 
model  is  still  superior  to  the  ECA  model. 

11.  Conclusion 

The  reliability  of  the  physics-based  SP  model  to  predict  the 
voltage  data  under  different  constant  current  charges  and  dis¬ 
charges  has  been  compared  to  the  semi-empirical  ECA  model.  The 
nonlinear  least  squares  approach  was  used  to  estimate  the  param¬ 
eters  of  the  models  by  minimizing  the  SSE  between  the  model 
prediction  and  the  data.  In  order  to  enhance  the  probability  of  find¬ 


ing  the  best  fit,  the  genetic  algorithm  was  also  implemented  to 
check  the  optimality  of  the  solution  obtained  by  the  trust  region 
method.  The  statistical  t- test  was  accomplished  to  determine  the 
significant  parameters  of  the  models.  For  the  SP  model,  the  cath¬ 
ode  diffusion  coefficient  did  not  influence  the  fitting  for  all  data 
sets.  The  anode  diffusion  coefficient  for  the  Quallion  data  and  the 
anode  reaction  rate  constant  for  the  USG  data  were  not  the  sig¬ 
nificant  parameters.  The  resistance  R  had  no  strong  effect  of  the 
fitting  for  the  ECA  model  for  all  data  sets.  Comparing  the  statis¬ 
tical  results  demonstrated  that  the  SP  model  fit  all  the  data  sets 
more  precisely  (less  SSE,  RMSE,  more  R 2,  AR2)  and  more  confi¬ 
dently  (tighter  prediction  bounds)  than  the  ECA  model.  On  other 
hand,  the  SP  model’s  computation  time  to  evaluate  the  voltages 
is  as  the  same  order  of  magnitude  as  the  ECA  computation  time. 
Thus,  the  SP  model  is  the  superior  candidate  model  rather  than 
the  ECA  model  and  it  can  be  applied  for  the  online  estimation  of 
the  SOC  and  the  state  of  the  health  of  the  battery  packs.  However, 
based  on  the  assumptions  of  the  SP  model  and  semi-empiricalness 
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of  the  ECA  model,  both  models  might  be  inaccurate  at  high  rates 
where  the  solution  phase  resistance  cannot  be  ignored  for  the 
SP  model,  and  more  electrical  components  for  the  ECA  would 
be  required.  The  other  limitation  of  the  present  models  is  the 
incapability  of  prediction  of  the  cell  capacity  fade.  Thus,  some 
extensions  (more  terms  and/or  parameters’  dependence  on  SOC 
or  time,  etc.)  to  the  current  models  are  essential  to  overcome  the 
restrictions. 
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